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ABSTRACT 

Classical T Tauri stars (CTTSs) are variable in different time-scales. One type of variability 
is possibly connected with the accretion of matter through the Rayleigh-Taylor instability that 
occurs at the interface between an accretion disc and a stellar magnetosphere. In this regime, 
matter accretes in several temporarily formed accretion streams or 'tongues' which appear 
in random locations, and produce stochastic photometric and line variability. We use the re- 
sults of global three-dimensional magnetohydrodynamic simulations of matter flows in both 
stable and unstable accretion regimes to calculate time-dependent hydrogen line profiles and 
study their variability behaviours. In the stable regime, some hydrogen lines (e.g. H/3, H7, 
HS, Pa/3 and Br7) show a redshifted absorption component only during a fraction of a stellar 
rotation period, and its occurrence is periodic. However, in the unstable regime, the redshifted 
absorption component is present rather persistently during a whole stellar rotation cycle, and 
its strength varies non-periodically. In the stable regime, an ordered accretion funnel stream 
passes across the line of sight to an observer only once per stellar rotation period while in the 
unstable regime, several accreting streams/tongues, which are formed randomly, pass across 
the line of sight to an observer. The latter results in the quasi-stationarity appearance of the 
redshifted absorption despite the strongly unstable nature of the accretion. In the unstable 
regime, multiple hot spots form on the surface of the star, producing the stochastic light curve 
with several peaks per rotation period. This study suggests a CTTS that exhibits a stochas- 
tic light curve and a stochastic line variability, with a rather persistent redshifted absorption 
component, may be accreting in the unstable accretion regime. 

Key words: line: profiles - radiative transfer - stars: variable: general - plasmas - magnetic 
fields - (magnetohydrodynamics) MHD 



1 INTRODUCTION 

Classical T Tauri stars (CTTSs) show variability in their light 
curves on time-scales from seconds up to decades (e.g. |Herbst et| 
|aL|2002l |Rucinski et"aTl|2008) . The origin of variability is likely 
different for different time-scales. The long-term variability may 
be connected with the viscous evolution of the disc (e.g. |Spruit| 
|& Taam||1993^ . The variability on the time-scales of a few stel- 
lar rotations may be connected with the inflation and closing of the 
field lines of the stellar magnetosph ere (e.g.|Aly & Kuijpers|1990| 
|Goodson, Winglee & Boehm||1997l|Lovelace et al.||1995), which 
is supported by the observations of AA Tau ( Bouvier et al.|2007] >. 
The variability on the time-scale of a stellar rotation are most likely 
connected with the rotation of a star, or with the obscuration of a 
stellar surface by a large-scale warped disc adjacent to the magne- 
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tosphere, which corotates with the star (e.g. |Bouvier et al.|[2003} 
|Alenca r et al. 2010; Roman ova et al.|2013| >. The variability on even 
smaller time-scales may be connected with the accretion of individ- 
ual turbulent cells formed in an accretion disc driven by magnetoro- 
tatinal instability (MRI) (e.g. |Hawley|2000"l|Stone et al.|2000|[R?>] 
|manova et al.|2012"| >. Another interesting possibility is that matter 
may accrete on to the star as a result of the magnetic Rayleigh- 
Taylor (RT) instability (e.g. |Arons & Lea| |l976; Sp ruit & Taam| 
[T99"3"l |Li & NarayanpOOl) . Global three-dimensional (3D) mag- 
netohydrodynamic (MHD) simulations (e.g. [Romanova, Kulkanu| 
|& Lovelace||2008| [Kuikarni & Romanova|2008| |2009) show that 
accretion proceeds through several unstable 'tongues', and the ex- 
pected time-scale of the variability induced by the instability is a 
few times smaller than the rotation period of a star. The variability 
on the time-scale less than a second is also expected if the inter- 
action of the funnel stream with the surface of the star leads to 
unstable radiative shocks which oscillate in a very short time-scale 
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Figure 1. Examples of the magnetospheric accretions in a stable (left panel) and an unstable (right panel) regimes. The background colours show the volume 
rendering of the density (in logarithmic scales and in arbitrary units). The sample magnetic field lines are shown as red lines. 



(e.g.|Koldoba et al. 2008 , S acco et al.|20T0) >. Finally, the variability 
can be also caused by magnetic flaring activities, which may occur 
on different time-scales (e.g. |Herbst et al.|1994"l|Wolk et al.|2005| 
|Feigelson et al.|2007| >. 

Here, we concentrate on the variability associated with the 
unstable accretion caused by the RT instability, and study spec- 
tral properties of stars accreting in this regime. The global 3D nu- 
merical simulations performed earlier (e.g. Romanova et al .|2008| 
Kulkarni & Romanova 2008 ) show that randomly-forming tongues 
produce temporary hot spots on the stellar surface with irregular 
shapes and positions, and that the corresponding light curves are 
very irregular. Such irregular light-curves are frequently observed 
in CTTSs (e.g. |Herbst et al.1[l994l |Rucinski et al^|200li) |Alencar| 
|et al.pOTO) . However, the irregular light curve can be also pro- 
duced, for example, by flaring activities. To distinguish the differ- 
ent mechanisms of forming irregular variability, we perform time- 
dependent modelling and analysis of emission line profiles of hy- 
drogen. 

To calculate the time-dependent line profiles from a CTTS 
accreting through the RT instability, we first calculate the matter 
flows in a global 3D MHD simulation with frequent writing of data. 
The fine time-slices of the MHD simulation data are then used as 



the input of the separate 3D radiative code TORUS (e.g. Harries 
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ton 2005 2006 ; Symington, Harries & Kurosawa 2005 ; Kurosawa 



Romanova & Harries 201 1 ; Kurosawa & Romanova 2012}. This al 
lows us to follow the time evolution of the line profiles that occurs 
as a consequence of the dynamical nature of the accretion flows 
in the 'unstable regime'. In earlier studies, we have used a simi- 
lar procedure (of 3D+3D modelling) to calculate line spectra from 
a modelled star accreting in a 'stable regime' (using only a single 
time slice of MHD simulations) in which the gas accretes in two 
ordered funnel streams (Kurosawa et al. 2008 1. More recently, we 
have performed a similar modelling for a star with realistic stellar 
parameters (i.e. V2129 Oph), and have found a good agreement be- 
tween the model and the time-series observed line profiles ( Alencar 
|et al.|2 012l. This example has shown that the combination of the 
3D MHD and 3D radiative transfer (3D+3D) models is a useful 



diagnostic tool for studying magnetospheric accretion processes. 
Examples of the magnetospheric accretion flows in both stable and 
unstable regimes are shown in Fig.[T] 

In this paper, we use a similar 3D+3D modelling method as 
in the earlier studies, but will focus on the CTTSs accreting in the 
'unstable regime'. A model accretion in a stable regime is also pre- 
sented as a comparison purpose. In the previous works, which fo- 
cused on the stable regime, we have fixed the magnetospheric ac- 
cretion at some moment of time (a single time slice of MHD sim- 
ulations) to calculate line profiles; however, in the current work, 
we will use multiple moments of time from MHD simulations to 
investigate variability of hydrogen lines caused by the dynamical 
changes in the accretion flows. In particular, we focus on the vari- 
ability phenomenon in the time-scale comparable with or less than 
a rotational period. 

In Section|2] we describe numerical simulations of stable and 
unstable accretion flows along with our numerical methods used in 
the MHD and radiative transfer models. The model results are pre- 
sented in Section [3] In Section [4] the model results are compared 
with observations. The discussion of the dependency of our models 
on various model parameters is presented in Section[5] Our conclu- 
sions are summarized in Section [6] Finally, some additional plots 
are given in Appendix [A] and in Appendix B in the online support- 
ing information. 



2 MODEL DESCRIPTION 

To investigate variable line spectra from CTTSs in stable and un- 
stable regimes of accretion, we perform numerical modelling in 
two steps. Firstly, we calculate the magnetospheric accretion flows 
using 3D MHD simulations and store the data at different stellar 
rotational phases. Secondly, we calculate hydrogen line profiles 
from the accretion flows, using the 3D radiative transfer code. For 
all the MHD and line profile models presented in this work, we 
adopt stellar parameters of a typical CTTS, i.e. its stellar radius 
R„ = 2.0 Rq and its mass M* = 0.8 M Q . Below, we briefly 
describe both models. 
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2.1 MHD models 

We calculate matter flows around a rotating star with a dipole mag- 
netic field with the magnetic axis tilted from the rotational axis by 
an angle O. The rotational axis of the star coincides with that of 
the accretion disc. A star is surrounded by a dense cold accretion 
disc and a hot low-density corona above and below the disc. Ini- 
tially, a pressure balance is present between the disc and corona, 
and the gravitational, centrifugal and gas pressure forces are in bal- 
ance in the whole simulation region (Ro manova et al.|2002[ l. The 
dipole field is strong enough to truncate the disc at a few stellar 
radii. To calculate matter flow, we solve a full set of magnetohydro- 
dynamic equations (in 3D), using a Godunov-type numerical code 
(see |Koldoba et al.||2002| |Romanova eTaLl[2Ufe] [20041 hereafter 
ROM04). The equations are written in the coordinate system that is 
corotating with the star. The magnetic field is split into the dipole 
component (Bd) which is fixed and the component induced by cur- 
rents in the simulation region (£>') which is calculated in the model. 
A viscosity term has been added to the code with a viscosity coef- 
ficient proportional to the a parameter (e.g. |Shakura & Sunyaev] 
11973b. The viscosity helps to form a quasi-stationary accretion in 



the disc. Our numerical code uses the 'cubed sphere' grid (Koldoba 
|et al.|20 02 l, and the Riemann solver used in the code is similar to 
that in Pow ell et al.| ( [T99 9). This code has been used for modelling 
the magnetospheric flows in stable and unstable regimes of accre- 
tion i n the past ([R omanova et al. 2003 , ROM04; Romanov a et al.| 
2008 ; Kulkami & Romanova 2008 1. It has been also used for mod- 
elling accretions to stars with a complex magnetic field (e.g. Long, 



Roma nova & Lov elace 2007 2008"l |Long et al.|2011[|Romanova et 
al.|2011b . 

In this work, we study the simulations performed for the sta- 
ble and unstable regimes of accretion (Fig. [TJ. Based on the previ- 
ous numerical simulations of magnetospheric accretions on to neu- 



tron stars by Kulkarni & Romanova (2008) (see also Kulkarni & 
Romanova 2005 1, we select two sets of model parameters that are 
known to produce a stable and a strongly unstable regimes. How- 
ever, we adjust some models parameters such that they become 
more suitable for a typical CTTS, and perform new MHD simu- 
lations. To keep a consistency, we adopt the same grid resolutions 
used in many of our earlier simula tions ( |Romanova et al. 2008 



Kulkarni & Romanova||2008| |2009| >, i.e. N r x N' z = 72 x 31 



grid points in each of six blocks of the cubed sphere grid, which 
approximately corresponds to N r x Ne x N$ = 72 x 62 x 124 
grid points in the spherical coordinates. 

The boundary conditions used here are similar to those in Ro- 
|manova et al.| (2003| > and ROM04. At the stellar surface, 'free' 
[<9(- ■ • )/dr — 0] boundary conditions to the density and pressure 
are applied. A star is treated as a perfect conductor so that the nor- 
mal component of the magnetic field does not vary in time. A 'free' 
condition is applied to the azimuthal component of the poloidal 
current: d(rB^,)/dr = such that the magnetic field lines have a 
'freedom' to bend near the stellar surface. In the reference frame 
that is corotating with the star, the flow velocity is adjusted to be 
parallel to the magnetic field B on the stellar surface (r = i?»), 
which corresponds to a frozen-in condition. The gas falls on to the 
surface of the star supersonically, and most of its kinetic energy 
is expected to be radiated away in the shock near the surface of 
the star (e.g. |Camenzind|1990|[K6nigl|1991[[Calvet & Gullbring| 
|1998| >. The evolution of the radiative shock above the surface of 
CTTSs has been studied in detail by |Koldoba et al.| ( |2008) . In this 
study, we assume all the kinetic energy of the flow is converted to 
thermal radiation (see ROM04 and Section [2/21 at the stellar sur- 
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Table 1. Basic model parameters used for the stable and unstable regimes 
of accretions. 



face. At the outer boundary (r = i? max ), free boundary conditions 
are applied for all variables. 

The key model parameters in our MHD simulations are: the 
stellar mass (M»), stellar radius (/?»), surface magnetic field on 
the equator (B aq ), corotation radius (r oor ), stellar rotation period 
(P*), tilt angle of dipole magnetic field with respect the disc axis 
(O) and a (viscosity) parameter. The parameter values adopted for 
the simulations in both stable and unstable regimes are summarized 
in Table Q] For more details on the simulation method and on the 
difference between the stable and unstable regimes, readers are re- 
ferred to |Romanova et al.|p003) , ROM04, [Romanova et"aT] ( 2008 ) 
and |Kulkarni & Romanova (20081. Although we have adopted a 
relatively large viscosity coefficient a = 0.1 for the model in the 
unstable regime (Table [jj, to be consistent with our earlier models 
( |Romanova et al.|2008||Kulkarni & Romanova|2008||2009> , our re- 
cent simulations have shown that the instability also develops with 
a smaller viscosity (e.g. a = 0.02-0.04). More detailed discussion 
on the dependency of our model on different model parameters is 
presented in Section]?] 



2.2 Radiative transfer models 

For the calculations of hydrogen emission line profiles from the 
matter flow in the MHD simulations (Section [2.1| l, we use the ra- 
diative transfe r code TORUS (e.g.|Harries|2000|[201 l[|Kurosawa et| 
|al.|2006[|201 l[[Kurosawa & R omanov a|2012| Lln particular, the nu- 
merical method used in the current work is essentially identical to 
that in Kurosawa et al. ( 201 1 ); hence, for more comprehensive de- 
scriptions of our method, readers are referred to the earlier papers. 
In the following, we briefly summarize some important aspects of 
our line profile models. 

The basic steps for computing the line variability are as fol- 
lows: (1) mapping the MHD simulation data on to the radiative 
transfer grid, (2) source function (S v ) calculations, and (3) ob- 
served line profile calculations as a function of rotational phase. In 
step (1), we use an adaptive mesh refinement (AMR) which allows 
for an accurate mapping of the original MHD simulation data onto 
the radiative transfer grid. The density and velocity values from 
the MHD simulations are mapped here, but the gas temperatures 
are assigned separately (see below). In step (2), we use a method 
similar to that of [ Klein & Castor ( 19781 (see also Rybicki & Ham-| 
|mer|1978[[Hartmann, Hewett & Calvet|1994b in which the Sobolev 
approximation (e.g. [Sobolev 1957, Castor 19701 is applied. The 
populations of the bound states of hydrogen are assumed to be in 
statistical equilibrium, and the continuum sources are the sum of 
radiations from the stellar photosphere and the hot spots formed by 
the funnel accretion streams merging on the stellar surface. For the 
photospheric contribution to the continuum flux, we adopt the ef- 
fective temperature of photosphere T p h = 4000 K and the surface 
gravity log g, = 3.5 (cgs), and use the model atmosphere of Ku- 
|rucz| ([T979). The sizes and shapes of the hot spots are determined 
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by the local energy flux on the stellar surface (see Section [372] for 
more detail). 

Our hydrogen model atom consists of 20 bound and a con- 
tinuum states. In step (3), the line profiles are computed using the 
source function computed in step (2). The observed flux at each 
frequency point in line profiles are computed using the cylindrical 
coordinate system with its symmetry axis pointing toward an ob- 
server. The viewing angles of the system (the central star and the 
surrounding gas) are adjusted according to the rotational phase of 
the star and the inclination angle of the system for each time-slice 
of the MHD simulations. For both stable and unstable cases, the 
line profiles are computed at 25 different phases per stellar rota- 
tion. We follow the time evolution of the line profiles for about 3 
stellar rotations; therefore, in total we compute 75 profiles for a 
given line transition. 

The gas temperatures in the accretion funnels from the MHD 
simulations are in general too high (T > 10 4 K) when directly ap- 
plied to the radiative transfer models, resulting in the lines being 
too strong. This is perhaps due to the 3D MHD simulations being 
performed for a pure adiabatic case (in which the equation for the 
entropy is solved). Solving the full energy equation would not solve 
the problem here because in both cases, the main issue is how the 
gas cools down at the inner disc and in the magnetospheric accre- 
tion flows. The process of cooling should involve cooling in a num- 
ber of spectral lines, and is a separate complex problem. Hence, 
we use the adiabatic approach for solving 3D MHD equations in 
which we solve the equation for the entropy. To control the gas 
temperature in the radiative transfer calculations and to produce the 
line strengths comparable to those seen in observations, we adopt 
the parameterized temperature structure used by Hartmann et al. 
l |1994} in which the gas temperatures are determined by assum- 
ing a volumetric heating rate proportional to r~ 3 , and by solving 
the energy balance of the radiative cooling rate ( Hart mann, Avrett] 
|& Edw ards 1982| and the heating rate. The same technique was 
also used b y e.g. | Muzerolle, Calvet & Hartmann ( 1998 20 0T) and 



|Lima et al. | < |20 1 0| >7 In our earlier work (e.g. [Kurosawa et al.|2008| ), 
the high temperatures from MHD simulations were 're-scaled' to 
lower values in order to match the line strengths seen in observa- 
tions. However, we found that the line profiles computed with the 
re-scaled temperature were very similar to those computed with the 
parameterized temperature structure from |Hartmann et al-H1994) 
(see Figure 7 in Kuro sawa et al.|20 08 1. Here, we adopt the latter as 
it is more physically motivated, but we expect to obtain similar line 
profiles even when we choose to use the re-scaled MHD simulation 
temperatures. 



3 RESULTS 

We perform two separate numerical MHD simulations for the mag- 
netospheric accretion as described in Section [2~T] with the param- 
eters shown in Table [T| which are known to produce a stable and 
an unstable magnetospheric accretion (e.g. Kulkarni & Romanova| 
2008 ). Typical configurations of matter flows in the stable and un- 
stable regimes are shown in Fig. [T] While matter accretes in two 
ordered funnel streams in a stable regime, it accretes in several tem- 
porarily formed tongues, which appear in random locations at the 
inner edge of the accretion disc, in the unstable regime. 

For the line and continuum variability calculations, we select 
the moments of time in the MHD simulations when the mass- 
accretion rates (M) of the system become quasi-stationary over 
a few rotation periods. Note that M for the unstable regime can 




Figure 2. The mass accretion rates (in 10 -8 Mq yr — 1 ) on to the central 
star plotted as a function of rotational phase for the stable (top) and unstable 
(bottom) cases. 



be still variable in a shorter time scale. The MHD simulation out- 
puts during 3 stellar rotation periods around these quasi-stationary 
phases are used in the variability calculations. The corresponding 
M as a function rotational phase (</>) for the stable and unstable 
cases are shown in Fig. [2] According to the figure, for the stable 
case, M rather steadily decreases by about 40 per cent during the 
3 rotation periods. On the other hand, for the unstable case, M in- 
creases slightly by about 20 per cent in 3 rotation periods, but it 
changes rather stochastically with smaller time-scales during those 
3 rotation periods. 



3.1 Line variability: persistent redshifted absorptions in the 
unstable regime 

The output from the MHD simulations are saved with the rota- 
tional phase interval of A(j> — 0.04, for 3 rotation periods. This 
corresponds to 75 total phase points at which line profiles are to 
be computed. This frequency is sufficient for 'catching' the main 
variability features in both stable and unstable regimes. For these 
moments of time, we calculate hydrogen lines profiles in optical 
and near-infrared i.e. Ha, H/3, H7, H<5, Pa/? and Br7 as a function 
of rotational phase, using the radiative transfer models as described 
in Section [2~2"l We set the gas temperatures in the accretion funnels 
for both the stable and unstable cases to be between ~ 6000 K and 
~ 7500 K. In all the variability calculations presented in this work, 
we adopt the intermediate inclination angle i — 60°. 

The subsets of the time-series line profile calculations for Ha, 
H/3, Hy, H8, Pa/3 and Bt - 7 are summarised in Figs. | A 1 1 and |A2| 
(in Appendix[A}, with A</> = 0.16, for the first rotation period. A 
quick inspection of the figures shows that with these particular sets 
of model parameters (see above and Table[T]l, we find no clear red- 
shifted absorption component caused by the accretion funnel flows 
in Ha profiles, and only a weak redshifted absorption component is 
seen in H/3 at some rotational phases. However, it is very notable in 
the higher Balmer lines (H7 and H<5), and also in the near-infrared 
lines (Pa/3 and Br7) partly because the peak fluxes in the emission 
part of the line profiles are relatively smaller than those of other 
lines (Ha and H/3). A similar tendency is also seen in some obser- 
vations (e.g. Edwards 1979; Appenzeller, letter & Jankovics 1986; 
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Figure 3. Magnetospheric accretion in the stable regime and the corresponding H<5 model profiles are shown at different rotational phases (indicated at the 
lower-left corner of each panel; with 0.08 interval). For a demonstration purpose, only a subset (only for about one rotation period of the star) of the MHD 
simulations and model profiles is shown. The volume rendering of the density shown in colour (in logarithmic scale) at each phase are projected toward an 
observer viewing the system with its inclination angle i = 60°. The model profiles are also computed at i = 60°. Line profiles models extended to the 
rotational phase (j> = 2.96 are presented in Fig. Bl (Appendix B in the online supporting information). 
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Figure 4. Same as in Fig.|3]but for the unstable regime. Line profiles models extended to the rotational phase i 
in the online supporting information). 
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: 2.96 are presented in Fig. B2 (Appendix B 



|Krautter, Appenzeller & Jankovics|199 0; Edwar dset al.|1994fr . As 
a demonstration for differentiating the line variability seen in the 
stable and unstable cases, we mainly use the evolution of the red- 
shifted absorption component in H<5 profiles. 



Fig. Ujj shows the 3D views of the matter flows (the volume 
rendering of the 3D density distributions as seen by an observer 

1 The animated version of this figure, extended to three rotation peri- 
ods, can be downloaded from |http : / /www .astro . Cornell ■ edu/| 
I -kur os awa/ research/ ctts_instab . html 
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with i — 60°) as a function of rotational phases, for the star ac- 
creting in the stable regime. The change in the mass-accretion rate 
during the first rotation is fairly smooth according to Fig. [2] and 
the geometry of the accretion funnels is also very stable. The cor- 
responding model line profiles of H<5 are also shown in the same 
figure. It shows a subset of the time-series calculations, i.e. only the 
phases between <j> — and cf> = 1.04 with the interval A(j> — 0.08. 
The time-series of model H<5 profiles during the whole 3 stellar ro- 
tations are given in Appendix B (Fig. Bl) in the online supporting 
information. Fig. [3] shows that the redshifted absorption compo- 
nent starts appearing around <j> ss 0, and its presence persists until 
<f> w 0.5. This phase span coincides with the time when the upper 
funnel accretion stream is located in front of the star, i.e., in the line 
of sight of the observer to the stellar surface. The redshifted absorp- 
tion occurs when the hotspot continuum is absorbed due to the line 
opacity in the accretion flow which is, for the most part, cooler than 
the hotspots and moving away from the observer. Between <f> w 0.5 
and 0.96, the redshifted absorption is not present because the lower 
accretion stream cannot intersect the line of sight of the observer 
to the hotspot. A similar behaviour was found in our earlier cal- 
culations with a stable accretion flow ( Kurosawa et al. 2008}. The 
appearance and disappearance of the redshifted absorption compo- 
nent is fairly periodic with their periods corresponding to the stellar 
rotation period (see also Fig. B 1 in the online supporting informa- 
tion.). 

A similar sequence of the 3D views of the matter flow in the 
unstable regime and the corresponding model line profiles for US 
are shown in Fig. [J] Again, the figure only shows a subset of the 
whole calculations, i.e. only for the phases between (f> — and 
(f> — 1.04 with the interval A(f> — 0.08. As before, the time-series 
of H<5 profiles computed for 3 stellar rotations are given in the sup- 
plemental figure in Appendix B (Fig. B2) in the online supporting 
information. The figure clearly shows that the geometry of the ac- 
cretion funnels and their evolution are very different from those of 
the stable case. The accretion no longer occurs in two streams, but 
rather occurs in the form of thin tongues of gas that penetrate the 
magnetosphere from the inner edge of the accretion disc (see also 
|Romanova et~aLl[2008l |Kulkarni & Romanova| |2008 ). The shape 
and the number (up to a several) of the tongues change within one 
stellar rotation period. The corresponding mass-accretion rate, dur- 
ing these rotational phases, changes rather stochastically as seen 
in Fig. [2] According to Figs.|4]and B2 (Appendix B in the online 
supporting information), the peak strength of the line also changes 
stochastically. Interestingly, in the unstable regime, the redshifted 
absorption component is present at almost all rotational phases. 
This is caused by the fact that there are a few to several accretion 
streams/tongues in the system at all times, and at least one of the 
accretion stream is almost always in the line of sight of the observer 
to the stellar surface. This behaviour is clearly different from that 
in the stable regime in which the redshifted absorption appears and 
disappears periodically with a stellar rotation. 

3.2 Photometric variability: stochastic variability in the 
unstable regime 

In both stable and unstable regimes, matter accretes on to the 
star and forms hot spots on its surface. In our earlier work 

2 The animated version of this figure, extended to three rotation peri- 
ods, can be downloaded from http :T/www . astro . Cornell . edu/| 
|~kurosawa/research/ ctts_instab . html 



(e.g. ROM04), the total kinetic energy of the gas accreting on to 
the stellar surface through the funnel streams are assumed to be ra- 
diated away isotropically as a blackbody. The effective temperature 
of a surface element in hot spots depends on the local energy flux 
at the surface element. The simulations show the hot spots are in- 
homogeneous with its highest energy flux and temperature in the 
innermost parts of spots. ROM04 have shown that the light curves 
in a stable regime are ordered/periodic, and have one or two max- 
ima per period depending on the tilt angle of the dipole and on the 
inclination angle of the system. On the other hand, in an unstable 
regime, the light curves are usually stochastic because the equato- 
rial accretion streams/tongues form stochastically and deposit the 
gas on to the stellar surface in a stochastic manner. This results in 
the stochastic formation of hot spots which produce the stochastic 
light curve (see also Romanova et al.|2 008 , Kulkarni & Romanova| 
2008). There is a wide range of possibilities when accretion is only 
weakly unstable in which both stable funnels and unstable tongues 
can coexist (Kulk arni & Ro manova 2009, Bach etti et aT1|2TJl0l >. 
However, in this paper, we focus on the case of a strongly unsta- 
ble regime as a demonstration. 

As mentioned earlier, our MHD simulations lack an imple- 
mentation of a proper cooling mechanism (e.g. radiative cool- 
ing), and consequently the adiabatically-heated gas in the funnel 
slows down gravitational acceleration of gas near the stellar sur- 
face. This results in the highest temperature of the hot spots to be 
only Tks ~ 4500K in the centre of the spot. This is not sufficient in 
explaining the ultraviolet radiation observed in many CTTSs (e.g. 
|Calvet & Gullbring|1998l |Gullbring et al.|2000| >. Here, we use a 
slightly different approach for the light curve calculation from hot 
spots. To overcome the low hot spot temperatures, we simply set 
the temperature to Th s = 8000K regardless of the local energy flux 
at the hot spots. Consequently, the temperature is uniform over the 
hot spots. The shapes and the locations of the hot spots are deter- 
mined by the energy fluxes computed at the inner boundary of the 
MHD simulation outputs. We apply a threshold energy-flux value 
such that the total hot spot coverage is about 2 per cent of the stellar 
surface (cf.fValenti, Basri & Johns|1993||Calvet & Gullbring|1998l 
|Gullbring et al.|2000||Valenti & Johns-Krull|2004| l. The threshold 
energy-flux values are fixed for each sequence of the line profile 
computations (in the stable and unstable regimes separately) so that 
the hot spot sizes can evolve with time when the local energy flux 
on the stellar surface changes. The continuum radiation from the 
hot spots is included in the previous line variability calculation in 
Section [3~T] No separate continuum radiative transfer model have 
been performed, as the continuum flux calculations are performed 
as a part of line profile calculations. The resulting variability of the 
continuum flux (light curve) near H<5 are shown in Fig[5]along with 
the hot spot distribution maps at four different rotational phases 
for both stable and unstable cases. The line wavelength of H<5 is 
4101 A which falls in a B-band filter; hence, the light curve shown 
here should be somewhat comparable to those from B-band pho- 
tometric observations. The same figure also shows the correspond- 
ing line equivalent widths (EWs) for H<5 as a function of rotational 
phase (cf. Figs. B 1 and B2 in Appendix B in the online supporting 
information). Similar variability curves for other hydrogen lines, 
i.e. Ha, H/3, H7, Pa/3 and Br7 are summarized in Appendix [X] 
(Fig.[A3). 

The surface maps for the stable regime in Fig.|5]show the pres- 
ence of banana-shaped hot spots on the stellar surface. For exam- 
ple, the hot spot on the upper hemisphere, which is created by the 
upper accretion stream (cf. Figs. [T] and [3}, is visible on the sur- 
face maps at the rotational phases (f> = 1.00 and 1.24. The spot 
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1.00 1.24 1.60 1.72 0.52 1.08 2.12 2.64 




Phase Phase 

Figure 5. The maps of hot spots, as seen by an observer at the inclination angle i = 60°, at four different rotational phases (top panels), the light-curves 
calculated at the wavelength of H<5 (middle panels; flux in arbitrary units) and the line equivalent widths (EWs) of model H<5 profiles (lower panels, in A) are 
shown for the stable (left panels) and unstable (right panels) regimes of accretions. The hot spots are assumed to be radiating as a blackbody with T = 8000 K, 
and are shown for representative moments of time which are marked as red dots in the light curves. 



on the lower hemisphere is visible, for example, at <fi = 1.60 and 
1.72. The light curve is periodic with its period corresponding to 
that of the stellar rotation. It is clearly seen that the maxima (at 
<f> w 0.25, 1.25,2.25) of the light curve occur when the upper 
hot spot faces the direction of the observer. The smaller secondary 
peaks at « 0.75, 1.75, 2.75 in the light curve are produced by 
the spot from the lower hemisphere. Interestingly, the line EW vari- 
ability curves almost mirrors the light curve, i.e. it shows a similar 
variability pattern but the maxima of the light curve corresponds 
to the minima of the EW values]^] This is caused by the combina- 
tions of the following: (1) for a given amount of line emission, the 
line EW decreases as the underlying continuum flux increases and 
(2) the amount of the redshifted absorption is largest when the ac- 
cretion funnel on the upper hemisphere is directly in between the 
stellar surface and the observer. Geometrically, this corresponds to 
a phase when the upper hot spot points toward the observer. The 
shapes of the spots are fairly constant during the three rotation pe- 
riods for which the light curves are evaluated. However, the peak of 
the light curves drop slightly in the later rotational phases since the 
mass-accretion rate decreases slightly during this time (see Fig. [2} 
and consequently the size of the spot decreases slightly. 

In contrast to the smooth and periodic light curve seen in the 
stable case, the light curve (at the H<5 wavelength) from the unsta- 
ble case is highly irregular, and it resembles many of the observed 
light curves from CTTSs (e.g., Herbst et al. 1994, see also Sec- 
tion |4.2} . The figure shows that there are typically several maxima 
per rotation period which correspond to several accretion tongues 
rotating around the star (Fig.[4}. The stellar surface maps show that 
the number of hot spots visible to the observer changes in time. For 
example, only one spot is visible at <f> — 0.52, and multiple spots 
are visible at <f> = 1.08 (two spots), 2.12 (two spots), 2.54 (3 spots 



3 In this paper, the sign of line equivalent widths (EWs) is opposite of that 
in a usual definition, i.e. a EW value is positive when the line is in emission. 



- the third one is barely visible but present on the upper right edge). 
The sizes of individual spots also change in time as this should cor- 
relate with the mass-accretion rate seen in Fig. [2] The same figure 
also shows that the variability of US EW is also stochastic, and no 
clear periodicity is found. Unlike in the stable case, the line EW 
for the unstable case does not clearly mirror the light curves, but 
they tend to correlate with each other. This is likely caused by the 
fact the mass-accretion rate and the density in the funnel flows are 
not constant in time for the unstable case while they are almost 
constant in the stable case. An increase in the mass-accretion rate 
would cause an increase not only in the continuum flux but also 
in the line emission because the density in the funnel flows should 
also increase for the higher mass-accretion rate. This effect would 
make the light curve and the line EW tend to correlate each other, 
and the mirroring effect seen in the stable case becomes weaken or 
disappears. The light curves computed at other line frequencies (at 
Ha, H/3, H7, Pa/3 and Br7) have also shown similar qualitative be- 
haviours, i.e. ordered curves in the stable regime of accretion, and 
stochastic curves in the unstable regime (see Fig. [A3}. In addition 
to the variable but rather persistent redshifted absorption compo- 
nent found in some lines (Figs. [4] |A2| and B2 - Appendix B in 
the online supporting information), the stochastic light curves and 
stochastic EW variability (Figs. [5] and |A3} are also key signatures 
of the CTTSs accreting in the unstable regime. 



4 COMPARISONS WITH OBSERVATIONS 
4.1 Line Variability 

Many spectral line variability observations of CTTSs, in different 



time-scales, have been carried out in the past (e.g.|Aiad et al. 



1984 



Edwa rds et al.|1994||To hns & Basri 1995a b; Gullbrin g et al.|1996 
Petrov et al.|1996t IChelli et al.|1997||Smith et all 19991 IPetrov etl 



|al.|199 9 , Oliveira et al. 2000 ; Petrov et al. 200T) |Alencar & Batalha| 



© 0000 RAS, MNRAS 000, 000-000 



Spectral variability ofCTTSs in unstable regime 9 



2002 ; Stempels & Piskunov 2002 ; Alencar et al. 2005 ; Kurosawa 



al. 2005 , Bouvier et al. 2007 , Donati et al. 2007 , Donati et al. 2008 



aet| 
)08J 



20T0l|2011a|bl | Alencar et al.|2012||Costigan et al.|2012||FaesTet| 
al.|2012 1. The line variability is often used to probe the geometry 



of magnetospheric accretion flows, and in many cases the accretion 
flows are found to be non-axisymmetric, e.g. SU Aur ( e.g. |Johns 



& Basri| 1995b) |Petrov et al.|1996| > and RW Aur A (e.g. |Petrov et 
al.|2001) . One of the key signatures of magnetospheric accretion 



is the presence of the redshifted absorption component in line pro- 
files. Such absorption components have been observed e.g. in some 



Balmer line s (e.g. |Aiad et al.||1984| [Appenzell er, Reitermann &| 
|Stahl| 1 988 ; Reip urth et al.|1996}, in n ear-infrared hydrogen lines 



(Pa/3 and Br7, e.g. 



Folha & Emerson 



lines (e.g. Beristain, Edwards & Kwan 



20011 and in some helium 



200T]|Edwards et al.|2 006). 



In some objects, observations show that the redshifted absorp 
tion component is rather persistently present in higher Balmer lines, 
e.g. in H7 and H5 in DR Tau and YY Ori (e.g. |Aiad et al.|1984| 
|Appenzeller et al.|[T988) |Edwards et al.||1994| >, in H<5 in DF Tau 
and RW Aur A (e.g. Edwards et al.|1994) . This contradicts with 
our stable accretion model (Fig. [3j in which gas flows on to the 
star occurs in two ordered funnel streams, and a redshifted absorp- 
tion is present only during a part of the whole rotational phase - 
when the funnel stream on the upper hemisphere is located in the 
line of sight of an observer to the stellar surface. However, a persis- 
tent redshifted absorption could be observed also in a stable accre- 
tion model when the system is viewed pole-on because the hot spot 
and the accretion funnel stream on one hemisphere is visible for an 
observer during a complete rotational phase (e.g. see Figure 5 in 
|Kurosawa et al.|2008) . 

No clear periodicity in the line variability is found in many 
CTTSs. For example, in the observation of TW Hya, Dona ti et al.| 
( |2011b) found that the variability of Ha and H/3 are not periodic 
(but see also Alencar & Batalha 20021, and suggested the cause of 
the variability is intrinsic (e.g. changes in the mass-accretion rates) 
rather than the stellar rotation. TW Hya also shows the redshifted 
absorption component in near infrared hydrogen lines (Pa/?, Pa7 
and Br7) in multiple occasions (e.g. |Edwards et al.|2006||Vacca &| 
|Sandell|2011) ; however, the frequency of the occurrence is unclear. 
In the studies o f DR Tau by | Alencar, Johns-Krull & Basri|p00T) 
and DF Tau by Johns-KrulT^sc Basri ( jl997| , no clear periodicity 
was found in their line variability observations. The non-periodic 
behaviour of the line variability is consistent with our unstable ac- 
cretion model (Figs. [4] and |5j. On the other hand, clear signs of 
periodicity are found in the line variability in some of CTTSs, 
e.g. BP Tau (e.g. Gullbrin get al.||1996[ [Donati et al.|[2008) and 



V2129 Oph (e.g. |Donati et al.|2007| |Donati et al.|2011a||Alencar 
et al.||2012] >. In the line variability study of V2129 Oph, Alencar 
et al. (2012) found that the redshifted absorption component in H/3 



was visible only during certain rotational phases. This is consistent 
with our stable accretion model (Fig. [3). 

4.2 Photometric Variability 

Observations show a variety of photometric variability in CTTSs 



(e.g. Bouvier et al.| 


1993] 


Gahm et al.||1993| |Herbst et al.||1994| 


Bouvier et al.|| 1995 


12003 


Stassun et al.||1999||Herbst, Maley & 


Williams 2000, Oudmaijer et al. 2001; Lamm et al. 2004, Alen- 


car et al.|2010|(. In some stars, clear periodic 


ight-curves are ob- 



served, while in others, the light curves appear stochastic (e.g., 
|Herbst et al~||1994[ |Alencar et"aT1|2010| >. Although the variabil- 
ity in X-ray are mostly stochastic, in some cases, signs of rota- 
tional modulation are found in soft X-ray emission in CTTSs (e.g. 




4540 4550 
HJD - 2,450,000 




2995 3000 3005 

JDH (2451545+) 



Figure 6. Top panel: the light-curve of TW Hya with the new broadband 
(located between V and R band) filter on MOST satellite (filled circle) 
overplotted with the V-band observation (open circle) of All Sky Auto- 
mated Survey (from Figure 4 in Rucinski et al. 20081. Bottom panel: A typ- 
ical stochastic light-curve (in 'white light' broadband) obtained by C0R0T 
satellite (from panel e of Fig. 1 in Alencar et al. 20101. Here, F and F0 de- 
note the flux and the maximum flux in the light curve, respectively. Credit: 
(upper panel) Rucinski et al. (2008), reproduced with permissions; (lower 
panel) [Aiencar et aTj^oTof Treproduced with permission © ESO. 



|Flaccomio et al.| 2005; Flaccomio, Micela & Sciortino 2012]). The 
well-ordered and periodic light curves are likely connected with 
rotation of hot spots resulting from ordered magnetospheric accre- 
tion streams (e.g., |Herbst et al!p 994), or from cold spots which 
represent regions of strong magnetic field. The exact origin of the 
stochastic variability is not well known, and there may be differ- 
ent causes. It was earlier suggested that stars with ordered light- 
curves are accreting in a stable regime while stars with stochastic 
light curves are in the unstable regime (e.g., Romanov a et al.|2008| 
|Kulkarni & Romanova|2008l[2009l l. 

Examples of periodic light curves can be found in e.g. Herbst 
[e7aI1p994) and |AlencareTaIJj2mO) . Using the COnvection Rota- 
tion and planetary Transits (C0R0T) satellite, [Alencar et al.H2010| l 
have shown that about 34 per cent of CTTSs in the young stellar 
cluster NGC 2264 display clearly periodic light curves with one or 
two maxima per stellar rotation (see panels a and b in their Fig. 1). 
Our model light curve from the stable accretion regime (Fig. [5} is 
very similar to the periodic light curves with two maxima per ro- 
tation. As in the observation, the smaller and larger maxima alter- 
nately appear. The periodic light curves with one maximum are also 
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consistent with a stable accretion observed with a lower inclination 
angle system (e.g. Figure 4 in Kurosawa et al. 20081. 

Examples of irregular light curves are shown in Fig. [6] The 
figure shows the light curve of TW Hya obtained with the Mi- 
crovariability & Oscillations of STars (MOST) satellite presented 
in |Rucinski et al.|p008|> and the CoR oT light curve of a CTTS in 
NGC 2264 presented in |Alencar et al.| |2010 ). Irregular light curves 
are found in about 39 per cent of the CTTSs sample in Alencar et al. 
( |2010} . |Rucinski et al.| ( |2008) finds no clear periodicity in the light 
curve for TW Hya in Fig. [6] However, they determined an approx- 
imate period of 3.7 d based on the light curve in an earlier epoch 
which showed some periodic behaviour. 

QlJsing 

the 3.7 d period, 

one can see that there are a few to several maxima per rotation in 
the light curve for TW Hya in Fig. [6] Interestingly, our model light 
curve for the unstable accretion regime (Fig. [5} also shows a few to 
several maxima per rotation. Similarly no clear period is found in 
the light curve of Ale ncar et al.| {2010) in Fig. [6} however, the pat- 
tern of the variability resembles that of our model from the unstable 
accretion regime (Fig. [5}. Further, the variability amplitudes of the 
irregular light curves in Fig.|6]are also comparable to the model in 
the unstable regime (Figs.[5land|A3l. 

In summary, we confirm that the periodic light curves found 
in the observations are consistent with a scenario in which the con- 
tinuum variability are caused by the rotation of two hot spots which 
are naturally formed in a stable accretion regime with a tilted dipole 
magnetosphere. Although not included in our model, the presence 
of stable cold spots on a stellar surface would also produce a peri- 
odic light curve. On the other hand, the irregular light curves found 
in observations are consistent with a scenario in which the contin- 
uum variability are caused by the rotation of stochastically form- 
ing hot spots which naturally occurs in the unstable regime. Yet 
the third type of the light curve (AA Tau-like), in which the vari- 
ability is partially caused by the occultation of a stellar surface by 
an accretion 'disc wall' or a warp near the disc truncation radius 



(e.g. Bouvier et al.|1999 Bertout 


2000 Bouvier et al.|2007||Alen-| 


|car et al.||2010| Romanova et al. 


2013), is not considered in this 



paper. 



5 DEPENDENCY ON MODEL PARAMETERS 

In general, the light curves and line EW variability behaviours are 
expected to be similar among the MHD models which clearly show 
stable and ordered accretion flows. The stable accretion funnels will 
always form smooth and periodic light curves and line variability. 
On the other hand, light curves and line variability behaviour for 
the unstable regimes will be slightly different, depending on the 
strength of the instability. In earlier models, (e.g. |Kulkarni & Ro-| 
manova 2008, Kulkarni & Romanova 2009), we have found that 
the number of unstable accretion filaments or 'tongues' present at 
a given moment of time decreases from several to just one or two 
as the instability changes from a strong regime to a weak regime. 
While the number of tongues may change for different model pa- 
rameters, the occurrences of the unstable tongues are still non- 
periodic. Therefore, the corresponding light curves and line vari- 
ability will also remain non-periodic, similar to those seen in the 
unstable model in this paper (Fig. [5}. However, the frequency of 
the peaks in the light curves and line EW curves per stellar rotation 
will decrease when the instability becomes weaker. Such variations 



This may be a quasi-period or drifting period (see also Siwak et al. 201 1 



are expected during the transition from a strongly unstable regime 
to a stable regime. The strength of the instability depends on some 
MHD model parameters such as P,, r cor , 6 and a in Table [TJ In 
the following we briefly discuss the dependency of our model on 
these key model parameters that can influence the stability of the 
accretion flows. 

Dependency on P* and r cor . These are essentially the same 
model parameters since the corotation radius (r cor ) is determined 
by the stellar rotation period P, for a given stellar mass (M*), 

1. e. r cor = (GAf*) 1/3 (P»/27r) 2/3 . The unstable flows occur more 
easily when the effective gravity (including the effect of the cen- 
trifugal force) is stronger at the inner rim of the accretion disc (e.g. 
|Spruit et al.|1995| l. A slower stellar rotation or equivalently a larger 
r cor will increase the effective gravity; hence, it tends to enhance 
the instability. The effective gravity becomes negative at the disc- 
magnetosphere boundary when r cor becomes larger than the trun- 
cation (magnetospheric) radius r m (where the matter stress in the 
disc becomes comparable with the magnetic stress in the magne- 
tosphere). Hence, in this simplified approach, accretion flows are 
expected to become unstable when r m < r cor (e.g. Spruit & Taam 
1 1 993[ >. However, the onset of instability depends on a few other fac- 
tors, such as the level of differential rotation in the inner disc, and 
also on the gradient of the ratio between the disc surface density 
and magnetic flux (e.g. Kaisig, Tajima & Lovelace 1992; Lubow & 
|Spruit 1995, Spr uit et al.|1995[ >. Hence, the condition, r m < r cor , 
is necessary for the onset of instability, but it may not be sufficient. 
In general, we find the instability is strong in our MHD simulations 
when stars are rotating relatively slowly i.e. r m < r cor . Conse- 
quently, we have chosen a larger P„ for the unstable case and a 
smaller P» for the stable case in this paper. 

Dependency on O. The instability can occur at different mis- 
alignment angles of the dipole O. However, at a large O (when the 
dipole is strongly tilted), the accretion flow becomes more stable 
because the magnetic poles are closer to the inner rim of the disc. 
The potential barrier in the vertical direction that the gas at the inner 
rim of the disc has to overcome to form the stable funnel accretion 
flow is reduced at a high O angle ( Kulkarni & Rom anova|2009| >. 
On the other hand, when O becomes small, the potential barrier to 
overcome becomes higher, and the stable accretion funnels cannot 
be formed so easily. In this case, the gas accumulates at the inner 
rim of the disc, and the accretion on to the star proceeds through 
instability 'tongues' which penetrate through the magnetosphere. 
In our MHD simulations, we find that a strong instability occurs 
when O < 25° (Kul karni & Rom anova 2009). Therefore, in this 
paper, we have chosen the large tilt angle Q = 30° to demonstrate 
the accretion in the stable case, and the small tilt angle O — 5° in 
the unstable case. 

Dependency on a. We use the viscosity parameter (a) to con- 
trol the global mass-accretion rate. The initial density conditions 
for the disc are common in all of our models, including the ones 
presented in the past (e.g. Kulkarni & Romanova 2008 1. The higher 
the value of a, the higher the mass-accretion rate becomes. At a 
higher accretion rate, the inner edge of the disc moves closer to the 
star. Consequently the ratio r m /r cor becomes smaller, and the con- 
dition becomes more favourable for the instability. The compres- 
sion of gas near the magnetosphere (the enhanced surface density 
per unit magnetic flux) may also play a role. The possible artifacts 
of the a viscosity on the formation of instability have been studied 
in a few experiments. For example, we start from a stable case with 
a small viscosity (a = 0.02), which is used for the stable case as in 
Table[TJ and we increase the initial density in the disc by a factor of 

2. In this case, we have observed that the accretion becomes unsta- 
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ble ( jRomanova et al.|2008| l. In another experiment, we start from 
a relatively small viscosity (a — 0.04), and increase the period of 
the star P, . In this case, we find the accretion becomes strongly un- 
stable (e.g. see Fig. 10 in Kulkarni & Romanova 2009). Moreover, 
in our recent simulations with the grid resolutions twice of those in 
the simulations presented in this paper (i.e. N r x iV 2 = 160 x 61 2 ) 
but with a smaller viscosity (a = 0.02), we find the accretions flow 
also become unstable, if a star rotates relatively slowly (e.g. see 
Fig. 17 in |Romanova et aL]2 013 1. Based on these experiments, we 
conclude that a large a is not a necessary condition for the insta- 
bility. For the MHD simulations presented in this paper, we have 
used a — 0.1 in the unstable case and a = 0.02 in the stable case 
for the consistency with earlier simulations (e.g. Romanova et al. 
[20081 |Kulkami & Romanova|2009> . 



6 CONCLUSIONS 

We have investigated observational properties of CTTSs accreting 
in the stable and unstable regimes using (a) 3D MHD simulations 
to model matter flows, and (b) 3D radiative transfer models to cal- 
culate time-series hydrogen line profiles and continuum emissions. 
In modelling line profiles, we have introduced the parameterized 
temperature structures of Hart mann et al.| (11994) in the accretion 
flow (Section |2.2fr and a fixed hot spot temperature (Section [3.2[ > in 
order to match the predicted line strengths and the amplitudes of 
light curves more closely to observations. Our approach is to intro- 
duce minimal modifications to the MHD solutions to obtain reason- 
able line profiles and light curves. Although this is not completely 
consistent with the MHD solutions, we still retain the flow geom- 
etry, density and velocity information from the MHD simulations. 
To implement proper heating/cooling mechanisms in the flows and 
more realistic hot spots, we need to develop proper physical models 
for them, which we have not done yet. 

In this work, we have focused on the variability in the time 
scale comparable with or less than a rotational period. For the un- 
stable accretion model, we have considered the flows associated 
with the Rayleigh-Taylor instability that occurs at the interfaces 
between the accretion disc and stellar magnetosphere. In particu- 
lar, we have analyzed qualitative behaviours of the variability in 
the redshifted absorption component in H<5 and the continuum flux 
to find key observational signatures that can distinguish CTTSs ac- 
creting in the two different regimes (stable and unstable). Our main 
findings are summarized in the following: 

(1) In the stable regime, the emission lines vary smoothly 
(Fig. [3](, and their line EW show one or two peaks per stellar ro- 
tation period (Fig. BJl. In the unstable regime, the EW of spectral 
lines vary irregularly and more frequently — showing several max- 
ima per stellar rotation period (Fig. [5}. 

(2) In the stable regime, the redshifted absorption is ob- 
served during approximately one half of the rotation period (when 
one of two accretion streams is located in front of an observer), 
and is absent during another half of the period (Fig. |3J. In the 
unstable regime, the redshifted absorption is observed at most 
of the rotational phases due to the presence of several unstable 
tongues/accretion streams which frequently move into the line of 
sight of an observer to the stellar surface (Fig. |4j. 

(3) In the stable regime, the continuum emission from hot 
spots (the light-curve) varies smoothly, and it has two peaks per 
stellar rotation period (Fig. BJ. In the unstable regime, the light- 
curve shows irregular variability with several peaks per period cor- 
responding to several unstable tongues (Fig. [5}. 



We find that the redshifted absorption components are in gen- 
eral more visible in the higher Balmer series, e.g. in H<5 and H7 
(see Figs. |Al| and |A2fr and some near-infrared hydrogen lines such 
as Pa/3 and Br7 than in the lower Balmer lines (Ha and H/3). Al- 
though the line variability analysis presented here is mainly based 
on US, similar conclusions can be reached if we use the other red- 
shifted absorption sensitive lines mentioned above. In addition to 
these hydrogen lines, He I A 10830 may be also useful for probing 
the accretion stream geometry of CTTSs. as demonstrated by Fis- 
|cher et al.|p008> . 

The work presented here may help to differentiate possible 
mechanisms of variability. Since the irregular light curves found in 
many CTTSs can be connected with different physical processes 
(e.g. magnetic flares, an accretion from a turbulent disc, an ac- 
cretion through an instability), the additional information from a 
line profile variability analysis may help to disentangle the differ- 
ent processes. For example, magnetic flares are not relevant to an 
accretion process; hence, no correlation is expected between a light 
curve and a redshifted absorption (which characterizes the magne- 
tospheric flows). If a light curve is irregular and the redshifted ab- 
sorption is persistent and non-periodic, the system may be accret- 
ing through an instability, which can occur either from a laminar, or 
from a turbulent disc. On the other hand, if a light curve is irregular 
and the redshifted absorption is periodic (seen only during a part 
of a rotational phase), the system may be accreting through well- 
defined accretion funnels originating from a turbulent disc, which 
provides turbulent cells to the accretion funnels and consequently 
produces an irregular light curve. 

We have investigated only the case of a strong instability 
in which most of the matter accretes in the unstable equatorial 
tongues. However, in reality, many CTTSs may be in a moder- 
ately unstable accretion regime (e.g. Kulkarni & Romanova j2009| > 
in which both ordered funnel flows and chaotic tongues are present. 
In these cases, both a persistent redshifted absorption component 
associated with the unstable accretion tongues, and ordered vari- 
ability pattern associated with the stable accretion streams can 
be present simultaneously. However, if the two ordered accretion 
streams from the stable flow components dominate, one can expect 
only one redshifted absorption per stellar rotation period as seen in 
the current work (e.g . Fig. [3} . Interestingly, a recent study by |Ro-| 
|manova & K ulkarni (2009) has shown that the unstable flows can 
be more ordered and corotate with the inner disc when the size of 
magnetosphere is small and comparable with the radius of the star. 
They have also shown that the formation of unstable tongues can 
be driven by the density waves formed in the inner disc region and 
by the interaction of the inner disc with the tilted magnetic dipole 
field of the star (Rom anova et al.|2012fr . 

Although we have concentrated on the line variability associ- 
ated with the magnetospheric accretion, some line variability could 
be also attributed to a variable nature of the wind. The winds from 
the disc-magnetosphere boundary (e.g. |Shu et al.|1 994) are usually 
associated with episodic inflations and openings of the field lines 
connecting the magnetosphere and the disc on a time-scale of a few 
inner disc rotati ons (|Aly & Kuij pers 1990; Goodson et al.|[T997| 
|Romanova et al.||2009||Zanni & Ferreira||2012| l. If the line emis- 
sion from a wind (either stellar or disc) is not significant, a variable 
wind would cause a variability mainly in the blueshifted absorption 
component. However, the emission from the wind could contribute 
non-negligibly to a total line emission (e.g. |Kurosawa et al.|2006 



Lima et al.|2"oTol|Kwan & Fischer 201 1 , Kurosawa & Romanova 



>va 



2012). If this is so, it makes the interpretation of a line variability 



a harder problem. On the other hand, some near-infrared hydrogen 
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lines e.g. Pa/3 and Br7, which show relatively strong redshifted ab- 
sorption in our models (Fig |A2} , are usually not affected by the 
winds (e.g. Folha & Emerson 2001); therefore, the diagnostic tools 
presented in this paper are still valid. In addition, the main time- 
scale of variability associated with winds is expected to be longer 
than that of the magnetospheric accretion considered in this paper. 
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APPENDIX A: ADDITIONAL LINE AND CONTINUUM 
VARIABILITY MODELS 

In the main text of the paper, we have used the subsets of H<5 model 
profiles and the continuum flux near US to demonstrate the key 
observable differences between accretion in the stable and unstable 
regimes. In this section, we present (1) the time-series profiles of 
additional hydrogen lines, i.e. Ha, H/3, H7, Pa/? and Br7 in the 
rotational phase (<j)) between and 0.96 with the interval A<f> = 

0. 16 (Figs. |AT| and \A2\ and (2) the line EW variability curves for 
the same additional lines as in (1), but for three rotational phases 
(between </> = and 2.96) with the interval A<j> = 0.04 (Fig. [A3}. 
The variability curves of the continuum flux (light curves) at each 
line frequency are also shown in the same figure. 

In Fig. |A1| Ha line profiles for both stable and unstable 
regimes show relatively small variability. No redshifted absorption 
component is seen in both cases. This is mainly because Ha line 
emission occurs in a wider region in the accretion stream, not just 
near the base of the accretion stream (near the stellar surface), and 
the strong emission in the line wings due to line broadening ef- 
fects (e.g. the Stark broadening, Muzerolle et al. 2001 1. The model 
profiles for H/3 and H7 show a very weak redshifted absorption 
component (below the continuum level) in both stable and unstable 
regimes at some rotational phases. Since the absorption component 
in these lines are very weak, they are not suitable for studying the 
differences between the stable and unstable regimes. Although not 
shown here, the redshifted absorption in these lines becomes more 
visible if the accretion funnel gas temperature is decreased slightly. 

The redshifted absorption components are much more promi- 
nent in US (in optical), Pa/3 and Br7 (in near-infrared) as one can 
see in Fig. |A2| for both stable and unstable cases. For the stable ac- 
cretion case, the redshifted absorption component is visible in these 
lines during about a half of a rotational phase (between <j> ~ and 
~ 0.5). These phases approximately coincide with the phases when 
the upper accretion funnel stream is in front of the star (see Fig. [3}, 

1. e. when the upper accretion funnel is in the line of sight of the 
observer to the stellar surface. A weak double-peak feature visible 
in Pa/3 and Br7 in the latter half of the rotational phase (</> ~ 0.5 
and ~ 1.0) is mainly caused by the 'lack of emission' from the gas 
moving at zero project velocity toward an observer. This is due to 
the geometrical effect (two distinct accretion streams) and the fact 
that the emission in Pa/3 and Br7 occurs near the base of accretion 
stream where the flow velocity is large (~ 200 km s^ 1 ). At some 
phases, the gas at the inner edge of the accretion disc could also 
intersect the line of sight between the observer and the system, and 
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causes a small amount of absorption near the line centre since the 
intervening gas in the inner most part of the accretion disc is just 
rotating with a Keplerian velocity; hence, its projected velocity to- 
ward the observer would be around zero. This can also contribute 
to the double peaked appearances of Pa/3 and Br7. For the unsta- 
ble case, the redshifted absorption component in US, Pa/3 and Br7 
is present at almost all phases. This is due to the geometry of the 
accretion flows (see Fig.|4](. The instability causes a few to several 
accretion streams to be present in the system at all times. Hence, 
at least one of the accretion streams or tongues is almost always in 
the line of sight of the observer to the stellar surface, and producing 
the redshifted absorption component. 

Finally, in Fig. |A3| we summarize the light curves (monochro- 
matic) computed at the frequencies around Ha, H/3, H7, Pa/3 and 
Br7 along with their line EWs as a function of rotational phase (for 
three rotational phases) for both stable and unstable cases. Note 
that the plots for H<5 are not shown here because they have been 
already shown in Fig. [5] in the main text (Section |3. 2} . While the 
light curves are very periodic for the stable case, they are quite ir- 
regular/stochastic for the unstable case. In the stable case, the light 
curves clearly show two maxima per period (except for the case in 
Br7). These are caused by the rotation of the two hot spots on the 
stellar surface created by the two accretion streams (see Figs.|3]and 
[5}. The line EW variablity curves are also influenced by the pres- 
ence of hot spots since they contribute significantly to the total con- 
tinuum flux. However, in the light curve of Br7 for the stable case, 
the relative amplitudes are rather small and it is more difficult to see 
the maxima. The smaller amplitudes are seen because the relative 
flux contribution of the hot spot (assumed 8000 K and ~ 2 per cent 
surface coverage) to the total continuum flux is much smaller at 
the wavelength of Br7 (2.165 \im). The accretion streams on to the 
stellar surface and the hot spots are formed rather stochastically in 
the unstable cases, resulting in the irregular light curves. See the 
main text in Section [3~2| for a further discussion on the light curves 
and the temporal variability in line EWs. 
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Figure Al. The time-series line profile models (solid) of Ha, H/3 and H7 in the stable regime (left columns) and in the unstable regime (right columns) shown 
for different rotational phases (between and 0.96). The rotational phase increases from the top to bottom panel, and the corresponding phase values are 
indicated in the leftmost column. In each panel, the line profile is compared with the mean profile (dashed, averaged over 3 rotational phases). All the profiles 
are computed with the system inclination angle i = 60°. 
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Figure A3. The temporal evolutions of the continuum flux (top panels, in arbitrary units) and the line equivalent width (bottom panels, in A) of Ha, H/3, H7, 
Pa/3 and Bry in both stable (left columns) and unstable (right columns) regimes are shown for three rotational phases. The plots for H<5 are omitted here since 
they have already been shown in Fig.|5]in the main text. The continuum fluxes are measured in the vicinity of the corresponding lines. In all cases, the system 
inclination angle i = 60° is used. 
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APPENDIX B: ADDITIONAL FIGURES 
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Figure Bl. The temporal evolution of model H<5 profiles (solid) for the stable accretion regime, computed for three rotation periods and with the system 
inclination angle i = 60°. The profiles are computed every 0.04 of a whole phase. In total, 75 profiles are shown. The corresponding rotational phases are 
shown on the upper-right corners of each panel. The mean profile (dotted) from the 3 rotational phases is also shown in each panel, for a comparison. 
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Figure B2. Same as in Fig. Bl but for the unstable regime. 



